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We study a quantum phase transition which occurs in a system composed of two impurities (or 
quantum dots) each coupled to a different interacting (Luttinger-liquid) lead. While the impurities 
are coupled electrostatically, there is no tunneling between them. Using a mapping of this system 
onto a Kondo model, we show analytically that the system undergoes a Berezinskii-Kosterlitz- 
Thouless quantum phase transition as function of the Luttinger liquid parameter in the leads and the 
dot-lead interaction. The phase with low values of the Luttinger-liquid parameter is characterized 
by an abrupt switch of the population between the impurities as function of a common applied 
gate voltage. However, this behavior is hard to verify numerically since one would have to study 
extremely long systems. Interestingly though, at the transition the entanglement entropy drops from 
a finite value of ln(2) to zero. The drop becomes sharp for infinite systems. One can employ finite 
size scaling to extrapolate the transition point and the behavior in its vicinity from the behavior of 
the entanglement entropy in moderate size samples. We employ the density matrix renormalization 
group numerical procedure to calculate the entanglement entropy of systems with lead lengths of 
up to 480 sites. Using finite size scaling we extract the transition value and show it to be in good 
agreement with the analytical prediction. 

PACS numbers: 03.67.Mn,73.21.La,72.10.Fk,71.10.Pm 
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I. INTRODUCTION 

There has been a recent flurry of activity relating 
entanglement entropy (EE)i (known also as the von 
Neumann entropy in quantum physics, and related to 
the Shanon entropy in information theory 2 ^—, and the 
Bekenstcin-Hawking entropy in the framework of black 
holes^— ), to quantum phase transitions (QPTs)&2 in 
condensed mater— ii£. The notion of of EE for a many- 
body system in a pure state arises when one divides it 
into two distinct regions: A and B. The entanglement 
between the subsystems A and B is measured by the EE 
Sa/b related to pa or pb, the reduced density matrix of 
regions A or B, respectively. 

Specifically, using the Schmidt decomposition, one can 
express any many-body pure state of the entire system, 
|\&), as the sum of two orthonormal basis sets of regions 
A ({\<f>A,i)}) and B {{\<f>B,j)}), such that 



QDr 



|*) =Ziai\<f>A,i)®\<f>B,i), 



(1) 



with real a.i > obeying £ia 2 = 1. This basis is closely 
related to the cigenbasis of the reduced density operators 
Pa/b = Ttb/^I*), i.e., 



PA/B = ^i<^\(f>A/B,i)((f>A/B,i\- 



(2) 



A unique measure of entanglement between the two re- 
gions A and B, is the von Neuman entropy of the reduced 
density matrix: 




Sa/b = -S,-a?ln(a|), 



(3) 



FIG. 1: The system: two single state quantum dots (impuri- 
ties) each of which is coupled to a different ID lead. Electrons 
may hop from the left (right) lead to the left (right) dot, but 
not between the dots. Thus, the two dots are coupled only 
by an electrostatic interaction. 



equivalent to the Shannon entropy of the squared 
Schmidt coefficients a 2 . This measure is also called the 
EE. Evidently S A = S B - 

Much effort was devoted to establishing the connection 
between EE and QPTs in many-particle one-dinrcnsional 
(ID) systems. It has been shown that in the critical 
regime of these models there are deviations from the cel- 
ebrated area law^, which states that the EE Sa should 
depend only on the surface area between regions A and 
B, hence it is constant (i.e., independent of system size) 
for ID systems. This is indeed the case when the system 
has a finite correlation length (when it is gapecUi). For 
cases in which the correlation length is infinite (i.e., in 
the critical regime) a logarithmic correction appears in 
the EE with a universal prefactor. Using conformal field 
theory arguments the EE of an infinite one-dimensional 
Luttinger Liquid (LL)«a and various corresponding spin 
chains was calculated, and the universal prefactor was 
related to the central charge of the underlying confor- 
mal field theory 3 -! 4 ^. Corrections to the EE due to finite 
siz o 13 ' 14 and the presence of a defect 1 ^— were also con- 



sidcrcd. For example, a static impurity embedded in a 
spinless LL, no matter how weak it is, will result (for re- 
pulsive electron-electron interaction :] 21 ' 22 ) in effectively 
severing the sample at the impurity location at low ener- 
gies, leading to a vanishing ground state conductivity. It 
is therefore not unexpected that this is manifested in the 
behavior of the EE, which tends to vanish for an infinite 
LL with a static impurity^. Dynamical impurities, e.g., 
an impurity with a resonant state which may fluctuate 
between an occupied and a vacant configuration, may on 
the other hand lead to a different behavior. Understand- 
ing the effect of dynamical impurities on transporter— 
as well as on thermodynamic*^— properties of ID LLs, 
has recently garnered much interest. 

In the present paper we define three tasks. The first is 
to analytically demonstrate that the system depicted in 
Fig-IU two impurities (quantum dots, QDs), each coupled 
to a different external lead, which are intercoupled only 
via electrostatic interactions (no inter-impurity tunnel- 
ing), undergoes a Berezinskii-Kostcrlitz-Thouless (BKT) 
QPT&2 as function of the LL parameter and the dot- 
lead interaction. This transition is similar in nature to 
the QPT recently predicted in the case where a quantum 
point contact is coupled to one of the impurities^ 5 -. The 
QPT manifests itself in the behavior of the population of 
both impurities. One of the impurities (corresponding to 
a broad level) swaps its population with the second (nar- 
row) impurity as a function of a common external gate 
voltage applied on the dots. Similar switching in differ- 
ent contexts have been studied extensively^ 6 - - — . By map- 
ping this problem onto that of a multiflavored Coulomb 
gas we have recently demonstrated that if the leads are 
non-interacting (regular Fermi-Liquids) the population 
switching is steep but not abrupt^, in agreement with 
previous studies^—. Once the external leads are LLs, 
we will show that this switching becomes abrupt at a 
critical value of the electron interactions in the lead (i.e., 
at a critical value of the LL parameter) or between the 
lead and the impurity. 

The second goal is to understand the EE behavior in 
this system, especially its sensitivity to the QPT. As we 
shall see, the system we consider here may be mapped 
onto an effective Kondo model 6 ^. When the leads are 
nonintcracting or weakly interacting (the LL parame- 
ter is close to one) the system corresponds to a anti- 
ferromagnetic Kondo phase, and therefore the two dis- 
connected parts of the system are nevertheless entangled, 
leading to a finite EE. On the other hand, for stronger in- 
teraction (lower values of the LL parameter) the system 
is in a ferromagnetic Kondo phase, in which no entan- 
glement between the two parts exists (even though they 
are still correlated), and therefore the EE is zero. The 
crossover between the finite value and zero at the critical 
interaction corresponding to the QPT should be sharp 
for an infinite system. 

Finally, as a third goal, we use EE to numerically study 
the properties of this QPT. Density matrix renormaliza- 
tion group (DMRG) 67 ' 68 is the most convenient numer- 



ical method to calculate the ground state properties of 
this interacting system, since EE appears naturally in the 
procedure. One might wonder though why not calculate 
directly the population of the impurities and see when 
population switching becomes discontinuous? This turns 
out to be extremely difficult since, in order to distinguish 
between a real discontinuity and a sharp transition, one 
must exceed length scales of order of the inverse relevant 
Kondo temperatur e 45 ' 61 " —. As we shall see, the latter is 
quite small even for noninteracting leads, and is further 
suppressed by LL correlations (until it vanishes at the 
transition point). One possible solution to this problem 
is to use DMRG with soft boundary condition s] 69 ' 70 . Here 
we shall show that by combining EE and BKT finite size 
scaling it is possible to obtain a reliable estimate of the 
critical interactions in the lead for manageable lead sizes. 
This paper is organized as follows: In the next sec- 
tion we introduce the model. Then in Sec. IHII we detail 
the mapping of our model onto a Kondo model, and il- 
lustrate the appearance of a BKT QPT. An analytical 
calculation of the dependence of the transition point on 
the parameters of the system is given. In the following 
section (Sec. HV|) , we calculate the EE using the DMRG 
numerical method, employing three different finite size 
procedures to determine the critical value of the interac- 
tion in the lead. All the different procedures give similar 
estimates of the critical value, which are in good agree- 
ment with theory in the range of parameters for which 
the theory is valid. We conclude with a discussion of the 
results obtained in the previous sections (Sec. |V|) . 



II. MODEL 

Here we consider a situation in which two single-state 
QDs (impurities) are coupled each to a different ID wires 
in the LL regime (see Fig.[lJ. The two impurities are cou- 
pled only by an electrostatic interaction. The spin degree 
of freedom is ignored (experimentally this corresponds to 
the presence of a strong magnetic field or ferromagnetic 
dots and leads which polarize the electron spin). The 
system is described by the following Hamiltonian: 



H = H, 



lead 



H 



lead 



+ H S 



impurities: 



where the Hamiltonian H( clld of lead 
right, respectively) is given by 



(4) 
L, R for left, 



#, 



lead 



-ft 



N t -1 



3 = 1 



c\ jCtj+i + H.c. 



JVf-l 



(5) 



„t 



here c\ „• is the creation operator of an electron on the )th 

site of lead I (whose total length is Ng), n^j = cj -c^j is 
the population of such a site, tg (JJg) is the hopping ma- 
trix element (interaction) between the sites in the lead, 



and vg is the electronic density in the lead. At half- 
filling (yi = 1/2, which we will always use in the nu- 
merical calculations) the system is in the LL phase for 
—2 < Ut/tt < 2;±i with the LL interaction parameter gi 
and the velocity of excitations vg given b y 12 ' 72 : 



(6) 

(7) 



2f, 



2cos- 1 {-U e /(2t e )Y 

7iVi-[£V(2*<)] 2 



2 cos- 1 [C/,/(2t £ )] 
The impurities Hamiltonian is 



-^impurities 

+ U' L (n L , d - 



£Ln L ,d+e R n R4 +U LR (n L ,d - \) {n R , d - 

) (til,i - vl) + U' R (n R:d - |) (n Rt i - v R ) 

+ *l4cl,i + t' R d R c R ,i + H.c. (8) 



where d L , R -. is the creation operator of an electron on the 

left (right) impurity, and rig, d — d\dg is the population of 
the ^-th dot. Here £l(R) is the left (right) impurity en- 
ergy (which may be varied by applying an external gate 
voltage; we henceforth assume that the latter has equal 
effect on both impurities, si = e\ — V g ), and Ul R is the 
Coulomb coupling between the two impurities. Each im- 
purity is coupled to the corresponding lead by both a lo- 
cal electrostatic interaction of strength Ug, and a tunnel- 
ing term parametrized by a hopping matrix elements t' e . 
The latter gives rise to level broadenings Tg = 7r|^| 2 po, 
with po being the local density of states on the last site 
of the lead at the Fermi energy [which is equal to l/(-7rf^) 
at half filling] . 



III. MAPPING ONTO KONDO 

Our analytic calculations are based on the Andcrson- 
Yuval mapping to a Coulomb gas 3 ^ -41 ' 45 ' 73- — . In this 
approach one expands the partition function to all or- 
ders in the dot-lead tunneling matrix elements t' L and 
t' R , and evaluates the resulting correlation functions at 
t' e = 0. Thus, the partition function becomes a sum 
over all possible imaginary time histories of tunneling in 
and out of each level, and can be cast in the form of a 
grand-canonical partition function of classical particles 
(charges) representing these hopping events. Because of 
the dot-lead interaction, each such process involves Fermi 
edge singularity physic a 82 ' 83 , hence is associated with a 
Fermi edge singularity exponent k fes , which is defined 
through the long-time behavior of the correlation func- 
tion (4(t)q,i(t)cJ 1 (0)^(0)) - r- K FB S a t t' e = 0. In 
addition, the level width Ti should be replaced by an ef- 
fective value r FES , which includes the prefactors in this 
correlation function. The value of k fes plays a crucial 
role in the following. It has been studied in our previous 
works 3 ^— ; here we will repeat our main findings. 

When ge = 1 (the lead is noninteracting) , we have the 

usual resonant level model, for which k fes = (l — 8l) 



TABLE I: The Fermi edge singularity exponent and the effec- 
tive level width. See the text for further details. 

Non-interacting lead Bosonization General expression 
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T, 



TV 

F e cos(5 e cS ) 



■k \tt\ po cos(5e), where 5e = tan 1 {irpoU' l /2) 



FES 



andF 

is the phase shift experienced by the electrons near the 
Fermi energy in lead £ due to the corresponding dot-lead 
interaction 8 - 2 -. In the general situation (gg not necessarily 
equal to unity), standard bosonization treatment yields 
k fes = I 1 ~ 9eUf/(TTV()] 2 /g e , ve being the velocity of the 
bosonic phase excitations in lead £, and L FES = -k \t'(\ po- 
Taking the limit gg — I we see that, while within a 
fermionic description k fes is expressed in terms of the 
phase shifts 5g, turning to a bosonized framework these 
8t are replaced by their leading order dependence on U' t . 
This is due to the linearization of the spectrum. Hence, 
the values of k fes and r FES in any particular model are 
renormalizcd by irrelevant operators not appearing in the 
Luttinger model. Boundary conformal field theory argu- 
ments show that Kp ES is related to finite size corrections 
to the spectrum of the lead with different potentials at 
its ends, which may be evaluated by analytical or nu- 
merical means. The results indicate 3 ^— that in general 
«fes = C 1 _ 2 9eS e cS /n) 2 /g e and r FES = n |i^| 2 p cos(S e eS ) 
for some effective phase shift S^ s e [— ir/2, 7r/2] which re- 
duces to the usual phase shift Si when the lead is nonin- 
teracting. For our model of the lead, i.e., a tight-binding 
chain with nearest-neighbor interactions, one can employ 
the Bethe ansataZ 2 . to find^i 1 .: 



Si s = tan" 



m 



VWe 



u 



(9) 



This discussion is summarized in table HI 

The Coulomb gas partition function can then be writ- 
ten in a standard form 7 -^— . The imaginary time history 
of the system (which is a circle whose circumference is 
the inverse temperature 1/T) is divided into intervals in 
which the system is in one of four possible states of the 
two dots: a = 00, 10, 01, and 11, corresponding to both 
dots being empty, only the left dot being occupied, only 
the right dot being occupied, and both dots being pop- 
ulated, respectively. The state a has a dimensionless 
energy h a , measured in units of l/£, where £ is a short- 
time (high-frequency) cutoff. The intervals are separated 
by hopping events, which are the classical Coulomb gas 
particles, as mentioned above. Their minimal separation 
is limited to £, which is thus of the order of the inverse 
bandwidth, £ ~ l/t£- A transition from configuration 
a to configuration (3 (a / /)) is associated with a fu- 
gacity y a p = y^ a , and a two-component vector charge 
e Q/ g = — e/3 a (the two components correspond to the left 



and right lead, respectively), obeying the triangle rule 
e Q7 + e* 7/ g = e a p. Physically, the components of e a p rep- 
resent the effective change in the charge of each lead in 
the corresponding transition. Values of these parameters 
for the system discussed are summarized in Table [TT1 The 
partition function reads: 



TABLE II: Parameters appearing in the Coulomb gas expan- 
sion, Eqs. (|10p - (|ll|l . The values of the Fermi edge singularity 
exponents Kfes and the effective level widths F FES are sum- 
marized in TablefT] 
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(10) 



where N+ 1 = 1, so that rjy+i = T\ + 1/T. The classical 
Coulomb gas action is: 



S({n,Oi})- J2 e ' Q ' Q -+ 1 ' e '^^+ lln { s in[7rT(T J -T l )] j 



iV 



7 j "-«;+! 



Tj+l — Tj 
I 



(11) 



One can now write down a set of 15 renormalization 
group (RG) equations for the Coulomb-gas parameters, 
which are valid to second order in the fugacities y a p but 
are otherwise exac t 73 ' 76 " —: 



dy a /3 _2 — n a p _ 
dln£ ~ 



-y, 



"/3 + H 



j/ Q7 j/ 7/3 e 



{h a +hp)/2-h~, 



(12) 



dn af 3 _ \-^ 2 h a -h a _V„2 phu-h-, n ri „\ 



dh a 
dln£ 



^-Ei^'^ + iEi^^. ( 14 ) 



/3,7 



where K a/ 3 = \e a p\ and k^ 7 = n a p + K ai - Kp-,. 

We will now concentrate on the Coulomb-blockade val- 
ley, i.e., |e |,£o + U > T L ,T R , where e = (e L +e R )/2, 
and we assume a small level separation el — £_r (see 
below). In this regime only the singly-occupied states 
are important at low energies. The RG flow is thus 
divided into three stages: (i) In the first one, £ _1 3> 
max(|eo|,£o + Ulr), and hence all the four filling con- 
figurations of the dots must be treated on equal footing, 
(ii) Then, as one enters the regime min(|£ |, £o + Ulr) <C 
£ _1 <C max(|£o|,£o + Ulr), the state with higher energy 
among the unoccupied and doubly-occupied configura- 
tions becomes higher than the cutoff and is discarded; 
(iii) Finally, for £ _1 <C min(|£ |,£o + Ulr), only the 
singly-occupied states 10 and 01 are left. In this last 
stage what remains is a Coulomb gas of only a single type 
of transitions. It is thus equivalent to the one originally 
derived by Anderson and Yuval for the single-channel 
anisotropic Kondo model 7 - 3 -, indicating the equivalence of 
the two systems. Under this mapping the two states 10 
and 01 become, respectively, the up and down states of 
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1 FESS 
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1 FESS 



s£ 



eoo.io = 



eoo,oi 



k fes ,0 



0, \l KpBg 



eoi,n 
eio.oi = 



K FES ,0 



hoa = 
hw - £hi 

hoi — £r£ 

hu = 

(£l+£r + U L r)£, 



J/io,oi = 

J/oo,n = ebo.n = 



"FES' V FES 



( Y K FES' Y K FES J 



the spin. The main effect of the two first stages of the 
flow is to establish the fugacity of the 10 ^ 01 transition, 
which is akin to the spin flip part (J xy ) of the Kondo ex- 
change coupling (via virtual processes through the dou- 
bly occupied state 11 and the unoccupied state 00). In 
addition, these two first stages lead to renormalization 
of the corresponding Coulomb-gas charge (related to the 
J z part of the Kondo exchange) and the energy differ- 
ence between these states, which is analogous to a local 
magnetic field B z along the z axis applied on the Kondo 
spin. 

Comparing the two Coulomb gases we can extract the 
parameters of the equivalent Kondo model. To the lead- 
ing order in max(r FES )/ min(|e , ^|, £i + Ulr) we find: 



PoJz =1 



"FES ~ "-FES 
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FES 1 FES 



Q 2 «< (NO Q2 K *(h + u LR }0 



=L,R 



£t 



£e + Ulr 
(15) 



PO'Jxy 



2-\/r FES r FES 



£«.=+«&, O (l £ ol0 QkL+k* (i £ o + U L r\0 



' ^FEST^FES 



FES ' ^FES 



eo 



B z =e L 



FES 



P, 



£o + U L r 

(16) 



1 FES 



P 2K n BB (\e R \0 - P 2K n ES ({sR + U L r}0 



(17) 



where P^(x) = T(l - /j,/2)/x 1 ~ fi/2 with T(z) the gamma 
function, Qfj,(x) = (1 — p/2)P^{x), and all the parameters 
of the original model refer to their bare values. Without 
interactions (apart form the Coulomb coupling between 
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FIG. 2: (Color online) Renormalization group flow of the 
anisotropic Kondo model in the J xy —J z plane. In the isotropic 
ferromagnetic case (J z = — J xy , red straight line) and all 
the region beneath it (blue lines) parameters flow towards 
weak coupling, whereas in the isotropic antiferromagnetic case 
Jz = Jxy, green straight line) and all the surrounding region 
(green lines) parameters flow towards strong coupling. 



the dots Ulr) one has k fes = 1, and these expressions 
reduce to those obtained in previous studie d 45 ' 61- — : 



T TL 
PoJz = 



1 



1 



1 



1 



7T \e L + U LR \e L \J it \£r + U L r \e r \ 

(18) 

2ynr^ ( i i \ 
Po J - = ~^r- [t^tu^ + w\) ' (19) 

D Vl , sl + Ulr , T R e R + Ul R 

B z = e L - e R In : — : 1 In : — : , 

7T \el\ n \e R \ 

(20) 

Since Ulr is typically of the order of the bandwidth ~ t(, 
the Q a (x) and P a (x) functions only change the corre- 
sponding terms by factors of the order of unity with re- 
spect to the case k fes = 1. Hence, the main effect of 
having k fes ^ 1 thus comes through the first term on 
the right hand side of Eq. (|T5l) . 

It should also be noted that due to renormalization 
effects, B z depends not only on the separation between 
the levels El — er, but also on their average position Eq- 
In particular, for \el — er\ < \Tl — Fr\ the sign of B z 
will change as a common gate voltage is varied across 
the Coulomb blockade valley. Correspondingly, the sign 
of the average spin projection (S z ) will also change, i.e., 
the two levels will swap their populatio n 45 ' 61 " — . 

The famous Kondo renormalization group 
equation: 



.06.7:5 



d(paJxy) 
dln£, 

d(poJz) 
dln£ 



= {poJxy){poJz), 
= (PoJxy) , 



(21) 
(22) 



imply that (cf. Fig. ^ for J z > —\J X y\ (including 
the case Ui = and U' e — 0) the exchange couplings 



grow under RG flow, and we are in the strong-coupling 
(antiferromagnetic-like) Kondo phase, where at low ener- 
gies the Kondo spin is strongly-coupled (entangled) into 
a singlet with the environment. This behavior sets in at 
energies below the Kondo temperature, which is para- 
metrically smaller than any other energy scale in the 
problem, and, in the generic anisotropic case (J z ^$> J xy , 
where renormalization of J z is negligible) , scales with the 
exchange couplings as TkS, ~ (po\Jxy\) 1 ^ PoJ ^ ■ The aver- 
age Kondo magnetization will be a smooth function of B z 
(i.e., the population is a continuous function of the gate 
voltage) with the scale set by the Kondo temperature at 
low energies (T <C Tk). For J z < — \J xy \ (which may 
occur only if k FFjS > 1) wc will enter the weak-coupling 
(ferromagnetic-like) Kondo regime, in which J xy flows to 
zero under RG, and the Kondo impurity becomes effec- 
tively decoupled at low energies, so that the population 
will be a discontinuous function of the gate voltage at 
zero temperature. The transition between the strong- 
coupling and weak coupling regimes (at J z = — \ J xy \) is 



of the BKT type, 
zero as ln(T^) ~ 



and the Kondo temperature goes to 
—{Jz + \Jxy\)~ 1 ^ 2 when approaching 
the transition from the antiferromagnetic side. 

This behavior can be understood in more physical 
terms. For example, having gi ^ 1 would shift pj z by 
l—(l/gL+l/g R )/2. By the Kondo renormalization group 
equations this indicates that J xy processes involving tun- 
neling out of one of the leads and into the other are more 
(less) relevant for gt > 1 (<% < 1), due to the enhanced 
(suppressed) tunneling density of states to the endpoint 
of a LL with attractive (repulsive) interaction o 12 ' 21 ' 22 . In 
addition, having U' t ^ modifies J z . The term linear in 
U[ in the expression for k FFjS represents the Mahan exci- 
ton effec t 82 ' 83 , and gives a positive (negative) contribu- 
tion to J z when U' e > (U' e < 0): when dot £ is occupied, 
then electrons arc repelled from (attracted to) the end- 
point of lead £, thus enhancing (suppressing) tunneling 
by the Pauli principle, and vice-versa. The quadratic 
term in U' t in the expression for k ffs represents the con- 
tribution of the Anderson orthogonality catastroph o 82 ' 83 , 
which always suppresses transitions and hence causes a 
decrease in J z . 

In the following numerical calculations we will confine 
ourselves to the case Ul = Ur = U (so that i?l = g R = g) 



^FES 



V FES 



kfes- We 



and U' L = U' R = U', i.e.. 
will also take Nl = N R = N. For el = er population 
switching will occur at the point of particle-hole symme- 
try, el = Er = —Ulr/2- The smooth-abrupt transition 
point will then be: 
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^Q 2K ; ES (^L«e/2) (23) 



U LR 



$ r HZ™ p Q^ u ^m 



These formulas are reliable for max(Fp Eg ) ^C Ulr, where 
the critical k fes is close to unity, so that Q2kJ es (Ulr£,/2) 
is also almost equal to one, and Tp ES rj T^. In addition, 



in the vicinity of the population switching wc will get 



B z 



r(2 - kfes) 



4(r^ ES 



r FEs)( £ 



Vlr/2) 



(UlrZ/2) 1 -^™ ttUlr 

(24) 
With the above expressions for ft FES one may find the 
critical U* and/or U'* for any value of Yi and Ulr- One 
may also test if the critical value U* at U' = and the 
critical value of U'* at U — are compatible (i.e., result 
in the same k fes ) even if Ti are not small enough with 
respect to Ulr, so that the critical k fes itself cannot be 
calculated analytically (we still need to assume r FES ra 



IV. NUMERICAL CALCULATION OF 
ENTANGLEMENT ENTROPY 

As discussed above, there are two possible phases 
for the system. One corresponds to the Kondo anti- 
ferromagnetic phase, while the other to the ferromagnetic 
phase of the Kondo model^. 

In order to calculate the EE we shall cut the system 
depicted in Fig.Q]into two parts, L and R. The left region 
(L) includes the left lead and the corresponding single- 
level QD (impurity), while the right region (R) includes 
the right lead and impurity. This is the most natural way 
to divide the system into two equal parts and has the 
great advantage of being very natural in the context of 
numerical DMRG calculations. Particle transfer between 
the two regions is prohibited. Nevertheless, since the two 
regions are coupled electrostatically, they are correlated 
and possibly also entangled. In the anti-ferromagnetic 
Kondo phase the ground state of the system is a singlet, 
i.e., an equal superposition of two states (one in which 
the left impurity is approximately full while the right one 
is approximately empty, and vise- versa) , resulting in two 
Schmidt coefficients a\ = 1/2 and a\ = 1/2, leading to 
S = S l /r = -ln(l/2) = ln(2). On the other hand, for 
the ground state of the ferromagnetic Kondo phase, the 
system is described by a single product many-body state, 
resulting in S = 0. The transition will thus manifest 
itself in a change of the EE between these two values. In 
the following we will utilize this change to facilitate the 
identification of the critical transition point. 



A. Dot-Lead Interactions 

We shall first consider the simplest case for which the 
leads are a Fermi liquid (i.e., U = hence g = 1) 
and the only interactions are the interdot interactions 
Ulr and the dot (impurity)-lead interactions U'. From 
now on we take the value of the hopping matrix el- 
ement t in the leads as our unit of energy. We set 
t' L = 0.1^2 and t' R = O.Ol-y/2, resulting in T L = 0.02 
and Tr = 0.0002, and vary Ulr and U'. Hereafter 
in this paper we tune the gate voltage to el — £r = 




FIG. 3: (Color online) Numerical computation of the EE S as 
function of the dot-lead interactions U' for different values of 
Ulr, the Coulomb interactions between the dots. The leads 
are noninteracting (U = 0). The hopping matrix elements 
between each lead and the corresponding dot are t' L — O.Ia/2 
and t' R = 0.01-^/2, resulting in Y L = 0.02 and F R = 0.0002. 
In all cases each lead is of length N = 40 (i.e., 40 sites). The 
transition between values of the EE close to ln(2) (indicated 
by dashed line) and to zero is quite sharp, and may be used to 
determine the critical value of [/'*. Lower inset: The ground 
state energy E(U') as function of the dot-lead interaction. A 
cusp in the energy, indicative of a level crossing, is clearly seen 
at values of U' corresponding to the sharp drop in the EE. 
Upper inset: The critical dot-lead interaction U'* as function 
of the interdot interaction Ulr- The symbols correspond to 
the numerical data while the line to the theoretical predictions 
of Eq. (|25|) . Good agreement is observed for large values of 
Ulr- 



—Ulr/2 — 0.0001, i.e., slightly above the expected pop- 
ulation switching point at half filling. In all forthcom- 
ing calculations up to 320 target states were kept for 
the longer systems. As detailed above, the mapping onto 
the Kondo model essentially depends only on r FES , r FES , 
and kfeS; which, for the present case, may be simplified 
to k fes = [l-2tan- 1 (C/72)/ 7 r] 2 . As depicted in Fig.O 
varying U' for a given value of Ulr will result in a modest 
increase in the EE as U' becomes more negative, up to a 
point where a sharp drop in the EE occurs. The larger 
the interdot interaction Ulr is, the smaller the absolute 
value of U' for which this drop occurs. As can be seen in 
the inset, this drop in the EE occurs for the same value of 
U' for which the ground state energy shows a downturn 
typical of a crossing between two distinct ground-states, 
as expected from a QPT. 

According to the analysis in the previous section, the 
transition should occur at [see Eq. (|2"3l ] 



U'* = 2 tan I - 



kUlr + 8yrlrR 
ttUlr - 4(r L + r fl ; 



1/2- 



(25) 



Comparing this prediction depicted in the upper inset of 
Fig. El by the curve to the numerical results for the tran- 




FIG. 4: (Color online) Numerical computation of the EE S as 
function of the dot-lead interactions U' for short lead length 
N < 120. The leads are noninteracting (U = 0). Here Ulr = 
1.2, and the hopping matrix elements between each lead and 
the corresponding dot are t' L — QA\/2 and t' R = 0.03\/2, 
resulting in Fl = 0.32 and Tr = 0.0018. Again, the transition 
between values of the EE close to ln(2) (indicated by dashed 
line) and to zero is sharp, and its location is quite independent 
of the system size, resulting in U' c as —0.47. Left inset: A 
zoom in of the transition region for larger system sizes 100 < 
N < 480. Right inset: the transition point U' c as function of 
system size. No substantial dependence of the critical dot-lead 
interaction on the system size is apparent. Thus, the infinite 
length critical dot-lead interaction is U'* = —0.47 ± 0.001. 



FIG. 5: (Color online) The numerically computed EE S as 
function of the interactions in the lead U for different lengths 
N of the leads. The dot-lead interactions are absent ([/' = 
0). The hopping matrix elements between each lead and the 
corresponding dot are t' L = 0.4\/2 and t' R — 0.03\/2 (Tl = 
0.32 and Tr — 0.0018), and Ulr = 1.2, corresponding to 
a Kondo temperature of T K « 9 • 10" 5 for U = Q&2k& 
In contrast to the dot-lead interaction case (Figs. EHIJ, the 
crossover between values of the EE close to ln(2) (indicated 
by dashed line) and zero is rather broad. As N increases the 
values of S for small U increase, while they decrease for larger 
values of the interaction, resulting in crossings of the curves 
for different values of S. Inset: Unlike for second order phase 
transition in which all curves for different N cross at the same 
value, for the BKT transition the crossing point shifts with 
N. 



sition corresponding to the symbols, one can see that the 
results match for large values of Ulr, although deviations 
appear at smaller values of Ulr- This is not surprising, 
since in order to obtain the relation depicted in Eq. ([25)) 
we assumed that max(Ti)/ULR <C 1. 

The system size considered in Fig. [3] is rather small 
(N L = Nr = N = 40). Nevertheless, one would not 
expect strong hnite-size effects on the critical dot-lead 
interaction U'* , since the mechanism by which tunneling 
is blocked (the Mahan exciton) is local, i.e., attraction 
of an electron in the lead to the vicinity of the occupied 
impurity due to the attractive lead-impurity coupling, as 
explained in the previous section. Indeed, as can be seen 
in Fig. 21 the main influence of increasing N is to bring 
the EE for —U'< —U'* closer to its expected infinite size 
value of ln(2), while the estimation of the infinite length 
dot-lead critical interaction [/'* as well as the sharpness 
of the transition do not seem to dependend on N. It 
is also interesting to note that S increases with \U'\ for 
—U'< —U'* (a trend which also appears in Fig.E)). This 
behavior could be understood as a result of the suppres- 
sion of the effective Fp ES by the dot-lead interaction (cf. 
Table UJ , consistent with the general trend observed in 
Fig. El 



B. Interactions in the Lead 

What happens if we consider a case with no dot-lead in- 
teraction (i.e., U' = 0), but instead with nearest-neighbor 
interactions in the leads (U ^ 0)? Since, as previously 
discussed, the phase transition does not depend on the 
details of the different interactions in the system, but 
rather on their combined contribution to the Fermi edge 
singularity exponent K FE gj22r— one expects that know- 
ing the critical value of U'* one is able to predict the 
critical value of U*. Thus, using Table U and Eq. © 
one can deduce that the critical LL parameter g* in the 
absence of dot-lead interaction will be connected to the 
critical value of U'* in the absence of interactions in the 
lead (U = 0, g = 1) via g* = (1 - 2tan([/'*/2)/7r)- 2 , 
from which we could find U* employing Eq. ([6]). This 
relation should hold even for T l,r/Ulr ~ 1, for which 
Eq. (|2"5j) no longer holds. Thus, for the parameters used 
in Fig. |H for which U'* = -0.47, one expects g* rj 0.75, 
resulting in U* rs 1. 

Indeed this can be demonstrated by the calculation of 
the EE depicted in Fig. [51 where the finite-size behavior 
of the same system considered in Fig. 0] (i.e., Ulr = 1.2, 
t' L = 0.4a/2, and t' R = 0.03^2, resulting in T L = 0.32 and 
Tr = 0.0018) is presented. A transition of S from values 



close to ln(2) to values approaching zero as U increases 
is apparent. Nevertheless, the transition is much more 
gradual than for the dot-lead interaction (see Fig. |4j, 
which makes sense since the LL character of the lead will 
be fully developed only for large sample sizes. A cross- 
ing between the curves corresponding to different sizes is 
apparent. While for small values of U larger system sizes 
correspond to larger values of S, the opposite occurs for 
large values of U, which is a hallmark of finite-size scaling 
of phase transitions 8 ^. 

It is important to emphasize that since we are dealing 
with a BKT phase transition, the finite-size scaling differs 
from that of a traditional second order phase transition 8 ^. 
While for a second order transition all the curves are ex- 
pected to cross at the same point, for a BKT transition 
the crossing point will drift as function of size and only 
at the limit of N — > oo will the crossing point correspond 
to the critical value U* of the lead nearest-neighbor in- 
teraction. The drift in the crossing point can be seen 
clearly in the inset to Fig. [5j where the crossing between 
the N = 30 and N = 40 curves occurs at U c rj 1.5, while 
between TV = 400 and TV = 480 U c « 1.05. 

Unlike for finite size scaling for second order transi- 
tions, there is no consensus on the optimal method to 
extract the transition point U*. Therefore, we shall em- 
ploy three different methods to identify the critical in- 
teraction for which the transition occurs: (I) extrapo- 
lation of the crossing point at N — > oo in the spirit 
of the phenomcnological renormalization group (PRG) 
procedur o 84 ' 85 ; (II) extrapolation of the transition point 
using a scaling ansatz inspired by the homogeneity condi- 
tion method 8 ^; (III) identification of the transition point 
by an heuristic scaling function. All these procedures 
give a similar estimate of U* , which is in good agreement 
with the correspondence between U* and U'* . 

For the PRG inspired (for details of the PRG proce- 
dure see Ref. I&H ) extrapolation of the crossing point as 
JV-^oo, one defines U c (Ni,N2), which is the value of 
the crossing between the curve corresponding to length 
iVi and the curve corresponding to N 2 , where N± and N 2 
are two successive (or next to successive) length values. 
We then extrapolate the behavior of the crossing points 
by the formula 8 ^ 



U C (N U N 2 ) = U* + 



C 



N1+N0 



(26) 



where C is a constant, and the value of the crossing point 
at infinite length is assumed to correspond to the tran- 
sition point U*. As can be seen in Fig. HI Eq. ([26]) 
works well for the entire range leading to U* = 1.02. 
When one examines more carefully the large N region 
(see the inset) the estimation of the critical point shifts 
a bit to U* = 0.96. Both values arc consistent with 
the expectations based on the dot-lead interaction data 
(U* = 1±0.1). 

The second method is inspired by the homogeneity 
condition method proposed and described in detail in 
Ref. |86|. In this method a function b(U, N) is constructed 
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FIG. 6: Obtaining the critical nearest-neighbor interaction in 
the lead U* by way of the extrapolation Eq. (|26l) . The symbols 
correspond to the crossing point between the curves in Fig. [5] 
for given values of Ni and JV2. Here we consider (Ni,N2) pairs 
that correspond to successive and next-to-successive values of 
N appearing in Fig. [S] The line represents a fit to Eq. (|26l) . 
Inset: A fit to Eq. (|26|) performed only for the large N values. 




FIG. 7: (Color online) B(U, N) [defined in Eq. (27J] as func- 
tion of the lead nearest-neighbor interaction U, for different 
system lengths N. The symbols correspond to the numerical 
results while the continuous curves to a fit using a 6th degree 
polynomial and the dashed curve to a fit using a 3rd degree 
polynomial in the crossing region 

. Inset: extrapolating the critical value U* for N — > 00. The 

symbols depict the values for which B(U, N) = 0, deduced 

from the curves appearing in the main figure. 



from the expectation value of the term in the Hamilto- 
nian driving the transition (in our case the interactions 
in the lead). Then the transition point is determined by 
the condition B(U*,N) = d N [N 3 d N b(U*,N)} = 0. Here 
we shall replace the expectation value b(U, N) by the EE 
S(U,N), and rewrite the function B(U,N) as a discrete 
differentiation: 

B(U, N) = N 3 S"{U, N) + [3/V 2 - (SN) 2 ]S'(U, TV), (27) 




FIG. 8: (Color online) Heuristic scaling [cf. Eq. J28J] of the 
different curves appearing in Fig. [5] The curves depict the 
scaling function S(U) = Sm ax (N)/S(U,N) for different sys- 
tem sizes N as function of U. Inset: Zoom into the transition 
region. Up to U ~ 0.95 all the curves coalesce. For larger val- 
ues the curves begin to separate. We identify this separation 
point as the critical U* . 



where S" (S") is the first (second) order discrete differ- 
entiation of S(U,N) {i.e., S'(U,N) = [S{U,N + SN) - 
S{U h N - SN)]/{2SN) and S"(U, N) = [S(U, N + 6N) + 
S(U, N-SN)- 2S(U, N)]/{SN) 2 }, and SN is the differ- 
entiation step. In Fig.[7]-B([/, N) is plotted for the larger 
system sizes (N = 160, 240, 320, 400) with SN = 80. 
Since for higher values of U there is strong scatter in 
the data (as indeed has been noticed in other applica- 
tions of the HMC, see cf. Ref. [861 ). we have interpo- 
lated B(U,N) over the whole range by a fit to a 6th 
degree polynomial (continuous curves), and in the vicin- 
ity of B(U, N) — by a fit to a cubic polynomial (dashed 
curve). The transition point for each size is determined 
by B(U, N) = 0, i.e., the point in which the curve crosses 
the x-axis. This point is found to be insensitive to the de- 
gree of the interpolation polynomial used. Extrapolating 
the critical value U* from the finite size values corre- 
sponding to B(U, N) = (see inset of Fig. [7]) results in 
U* = 0.94. 

The third method we shall use in order to determine 
the critical value of U* involves a heuristic scaling func- 
tion of the different curves appearing in Fig. [5] This 
scaling is based on the observation that the maximum 
entropy for all system sizes seems to appear at the same 
value of U « 0.55, and that the general form of S(U, N) 
around this value shifts in a rigid manner as function of 
N. Thus, we postulate that for U < U* one may col- 
lapse all the different curves onto a single curve by the 
following heuristic scaling function: 



S(U) = 



Smzx(N) 



(28) 



of the EE for a given length N of the system. As can be 
seen in Fig. [FJ applying this heuristic function collapses 
all the different curves on a single one up-to a certain 
value of U, while for larger values the curves diverge. 
This is somewhat similar to the situation one encoun- 
ters for the PRG scaling of the energy gap in the BKT 
transition^. Identifying the point for which the curves 
begin to diverge with the critical point (as is done in the 
PRG procedure) results in U* ~ 0.95, in agreement with 
the other methods we employed here for identifying the 
critical point. 



C. Both Dot-Lead and Intra-Lead Interactions 

Now we shall discuss the behavior for the case where 
both U' and U are nonzero. We will consider the case 
where both U' > and U > 0, which is the case ex- 
pected to be encountered in realistic experimental de- 
vices. We shall set the dot-lead tunneling matrix ele- 
ments t' L = 0.2V2 and t' R = 0.02V2, corresponding to 
T L = 0.08 and T R = 0.0008, while retaining U LR = 1.2. 
For these parameters ~nieix(Ti)/ULR <C 1 and we expect 
the determination of the critical point which depends on 
the parameter k ees given in Eq. (|23[) to hold. Setting 
V = 0.3 and using the relationship given in Eq. fl5J, the 
transition is expected to occur at U* obeying the follow- 
ing relation: 



U' = ^4-(C/*) 2 x 



tan 



cos 




(29) 



S(U, N) ' 
where S max (N) = S(U — 0.55, N) is the maximum value 



The right hand side, for a given value of T^, T#, and Ulr, 
depends only on U*. Thus by plotting the right side of 
Eq. (|2"9")l and determining where it crosses the value of 
U', it is possible to evaluate U*. For the parameters con- 
sidered here, the critical value corresponds to U* = 0.96. 
There is an additional crossing at U* = 1.95, but since 
this is very close to the LL-charge density wave phase 
transition point (at U = 2) it would be extremely hard 
to observe it numerically for reasonable system sizes. 

One may wonder why in the previous sub-section we 
did not consider also the case of max(r£)/C/L_R <C 1, 
which would facilitate a direct comparison between U* 
and theory, instead of discussing a case for which r^ 
are comparable to Ulr, leaving only the possibility to 
compare between the numerically computed U'* and U*l 
The reason is that for max(r^)/{7_L_R -C 1 and V = 0, U* 
tends to be close to zero, i.e., g w 1. For values of the 
LL close to the non-interacting case, one must go to very 
large systems to see the LL behavior developing^, which 
is beyond our current capabilities. 

Returning to the case at hand, it can be clearly seen 
in Fig. IH1 that as in the previous cases, a typical BKT 
crossing of the S as function of U for the different lengths 
is observed. Using the extrapolation of the crossing point 
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FIG. 9: (Color online) The numerically computed EE S as 
function of the interactions in the lead U, for a constant dot- 
lead interaction U' = 0.3 and different lengths N of the lead. 
The hopping matrix elements between each lead and the cor- 
responding dot are given by t' L — 0.2\/2 and t' R — 0.02\/2 
(T L = 0.08 and F R = 0.0008), while Ulr = 1.2, correspond- 
ing to a Kondo temperature of T K « 10" 13 at U = o^SLffi 
The transition between values of the EE close to ln(2) (indi- 
cated by dashed line) and to zero is steeper than in Fig. [5] 
although the typical BKT transition crossing are still clearly 
evident. Inset: An extrapolation of the critical U* in the 
same manner as in Fig. [6] The symbols depict the crossing 
point between the curves in the main figure corresponding to 
Ni and N2, where successive and next to successive lengths 
are taken. The line represents a fit to Eq. (|26[) . 



at N — > 00 following the PRG procedure, and fitting it 
to Eq. ([26]) (see inset), an extrapolated value of U* = 
0.98 is obtained, in good agreement with the analytical 
prediction. 



system corresponds to a ferromagnetic Kondo model, for 
which there is no entanglement between the sub-systems. 
Using this behavior of the EE, and the fact that it lends 
itself to straightforward calculation within the frame- 
work of numerical DMRG, we were able to use finite 
size scaling to identify the QPT between the smooth and 
abrupt switching phases. We have found that using the 
signature of the QPT on the EE behavior indeed gives 
an accurate method to study the transition properties. 
This reproduces the analytical results obtained by map- 
ping the the system onto a Kondo model. One may ex- 
pect that the EE could be used in a similar fashion to 

identify and study different QPTs related to dynamical 
impurities^>22^>29-ii£8 

Finally, it should be noted that according to Eq. (J23I) . 
in the absence of dot-lead interactions (If = 0), and for 
strong intcrdot interactions Ulr ^S> max(r^), the abrupt 
population switching will occur at values of the LL pa- 
rameter smaller than a critical value which is quite close 
to unity 



1-9* 



max(r£) 
Ulr 



(30) 



Thus, for large intcrdot interaction Ulr 3> maxT^, even 
a weak interaction in the leads will result in an abrupt 
population switching. Once dot-lead interactions are 
taken in account, the critical g* will be further decreased 
by 2U' /(irb). Since t is essentially the band width and 
the dot-lead interaction is parametrically much smaller, 
the abrupt population switching phase will appear for 
values of g « 1 even when repulsive dot-lcad-intcractions 
are included, and should thus be accessible experimen- 
tally. The population switching can be probed using a 
quantum point contact as a charge sensor^ 



V. DISCUSSION 

We have shown that in the presence of repulsive inter- 
action in the lead (U > 0, g < 1), the system depicted 
in Fig. [T]may show an abrupt population switching, i.e., 
an abrupt swap of the left dot-right dot population as 
function of the applied gate voltage. The nature of the 
population switching has a clear signature in the behav- 
ior of the EE. For smooth switching the system can be 
mapped on the antiferromagnetic Kondo model, resulting 
in a finite entanglement between the left and right sub- 
systems. On the other hand, for abrupt switching the 
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